Dynamical simulation of integrable and non-integrable models in the Heisenberg picture 
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The numerical simulation of quantum many-body dynamics is typically limited by the linear growth of en- 
tanglement with time. Recently numerical studies have shown, however, that for ID Bethe-integrable models 
the simulation of local operators in the Heisenberg picture can be efficient as the corresponding operator-space 
entanglement grows only logarithmically. Using the spin- 1/2 XX chain as generic example of an integrable 
model that can be mapped to free fermions, we here provide a simple explanation for this. We show furthermore 
that the same reduction of complexity applies to operators that have a high-temperature auto correlation function 
which decays slower than exponential, i.e., with a power law. This is amongst others the case for models where 
the Blombergen-De Gennes conjecture of high-temperature diffusive dynamics holds. Thus efficient simulabil- 
ity may already be implied by a single conservation law (like that of total magnetization), as we will illustrate 
numerically for the spin-1 XXZ model. 

PACS numbers: 02.70.-c, 67.80.dk, 75.10.Pq 



White's Density Matrix Renormalization Group (DMRG) 
IH and it's more recent generalizations to time evolution us- 
ing the Time Evolving Block Decimation (TEBD) J2l or t- 
DMRG [3] algorithms are indispensable tools in the numer- 
ical simulation of one-dimensional quantum many body sys- 
tems. They permit high-accuracy calculations, provided that 
the entanglement between any two complementary partitions 
remains small. For finite-range interactions this is the case 
for the ground state J3l. However in real time evolution the 
entanglement often grows linear in time, limited only by the 
Lieb-Robinson upper bound J5l|6[]. E.g. for the spin-i XY 
chain the evolution of the entanglement entropy was investi- 
gated in [Q] showing explicitly the linear growth in time. 

However the evolved state contains a lot of information 
which is of little interest. Experimental measurements as well 
as theories are almost solely concerned with few particle prop- 
erties, i.e., quantities that can be expressed in terms of expec- 
tation values of only a small number of elementary operators. 
This suggests to go to the Heisenberg picture (HP) instead, 
and to simulate the dynamics of these operators. Prosen et 
al. JS 0] where the first to pursue this approach. They ob- 
served an exponential speed up in numerical simulations of 
local operators for integrable systems. So far there is however 
no general understanding of why this is the case and whether 
or not integrability is crucial. In the present paper we provide 
an explanation of the speed-up for integrable models that can 
be mapped to free fermions. We also argue that integrability 
is not necessary and that the existence of a conservation law 
may suffice for the efficient simulation of the dynamics of lo- 
cal operators that constitute the conserved quantity. We will 
discuss the spin-i and spin-1 XXZ models as specific exam- 
ples supporting and illustrating our arguments. 

In order to do HP simulations using e.g. the TEBD scheme, 
the operator <D(t) at time t is expressed in terms of a ma- 
trix product operator (MPO). For typical observables this is 
straightforward for the initial time t = 0. Time evolution 



is then calculated by updating the matrices according to the 
Heisenberg equation of motion using a Trotter decomposi- 
tion. Efficient simulation requires that the matrix dimension 
of the MPO's (called bond dimension) is limited to a max- 
imum value x- This means that only the \ largest Schmidt 
values in the Hilbert space of operators are kept, correspond- 
ing to a small operator-space (OS) entanglement between any 
two complementary partitions of the lattice. To quantify the 
entanglement of an operator 0(t), which after proper normal- 
ization can be viewed as state vector in OS, we use the OS 
Renyi entropies (OSRE): 

log 2 Tr k 



S a — 



1 - a 



> Sr, (3 > a > 



(1) 



Here k is the corresponding reduced density matrix in OS re- 
sulting from tracing out the left or right partition at a given 
bond. In the limit a 1, S a is the well known von Neu- 
mann entropy, which is a good measure of bi-partite entan- 
glement. For a ~ > 0, S a gives the dimension of the Hilbert 
space. Clearly for an MPO of bond dimension x> me max- 
imum for all Renyi entropies is log 2 X- Although it is not 
yet fully established when a quantum state or an operator is 
faithfully represented by a matrix product with finite bond di- 
mension, one can employ the results of Schuch et al. ifioll 
to show that efficient simulation is impossible if the Renyi 
entropies with a > 1 scale faster than logarithmically with 
time. If S a> i grows linearly in time, we must expect that 
the computational cost required to reach a certain accuracy, 
which is polynomial in x< will grow exponentially with time 
(note that this is not necessarily true for S a with a < 1 (10]). 
In fact for the time evolution of typical state vectors in the 
Schrodinger picture this is very often the case ll 111 . On the 
other hand an at most logarithmic growth of S a> i is a nec- 
essary condition for an efficient simulability. Although not 
sufficient, it also gives good indication when such a simula- 
tion is possible. In the following we will discuss the time 
evolution of the OSRE for a generic model, the XXZ chain, 
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Aaj&j +1 ), where a x 



denote the Pauli matrices in the spin-i case and the spin-1 
matrices (eigenvalues —1, 0, 1) in the spin-1 case respectively. 
The spin-i case is integrable for any value of the anisotropy 
A. For the special case of A = (spin-i XX model) this 
model can be mapped to free fermions. 

integrable models equivalent to free fermions: Let us con- 
sider the spin-i XXZ model as a generic example of a ID 
integrable model. We have calculated the time evolution of 
the OSRE 52 for different types of simple operators using the 
TEBD scheme with open boundary conditions and a fourth 
order Trotter decomposition 11211 . The restriction to open 
boundary conditions is not an issue for local operators as long 
as the time is shorter than the propagation time to reach the 
boundaries [5]. Although not shown the OS von-Neumann 
entropy Si has the same scaling behavior. One clearly notices 
that the OSRE of all operators scales at most logarithmically 
in time, an observation made already by Prosen et al. for other 
integrable models In the special case of A = the en- 

tropy even saturates at a finite value for some operators like 
a z or products at a small number of different lattice sites. 
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FIG. 1: (Color online) OSRE dynamics for the 40 site spin-± XXZ 
model for a split in the center. The legend gives initial operator and 
anisotropy in the order in which the arrow cuts the graphs. Dashed 
lines mark infinite index operators (see text). \ — 1000 is used in all 
cases and the numerical error is negligible on the time scale shown. 

In the following we will provide an explanation of the en- 
tropy scaling for the case of the XX-model, i.e., for A = 0, 
which can be mapped to free fermions. This will be done by 
reexpressing the XXZ model in terms of Majorana-fermion 
operators 11411 . which turns out to be more convenient than 
the more common Wigner-Jordan transformation: w 2 j-i = 




U<j"ij a j' ™V = [UtejOfjOj- The Majorana op- 
erators are Hermitian and fulfill anti-commutation relations 
{uij, wi} = 25 ji . The three types of interactions in the XXZ- 
model can be reexpressed as 



j+1 ■= -W 2 jW2(j+l)-l, 

j+1 = iw 2 j-iw 2 ( j+1 ), 

1+1 = -W2j-lW 2 jW 2 U+l)-lW2(j +1 ) 



(2) 



A complete basis in the OS is given by P a = JT . w 2 j-i w 2 j 3 , 
where a = (ai,a 2 , ■ ■ ■ ) and {a;} <G {0, 1}^. We can now 
define adjoint-fermion annihilators and creators via a,j\P a ) = 
aj\WjP a ), a]\P a ) = (1 - aj)\WjP a ), with {a^a]} = 
Sji. Associating the adjoint vacuum |Pq) with the unity 
operator 1, i.e., |1) = |Po), we can express all opera- 
tors in terms of adjoint-fermion excitations lfl3ll : \P a ) = 



YlMi-i 



Mapping the Heisenberg equa- 
tion then gives a Schrodinger like equation for the evolution 
in OS, 



d - 
l dt P * 



P a ,H 



=-n\p a 



(3) 

with a "super"-Hamiltonian T-L. Explicitly calculating 
the terms in the commutator for the XX-model via 

\[P a ,a]a] +1 ]) = 2i(a|j02(j+i)-i - h.a.)|P Q ), and 

\[Pa,d-jO-j +1 }) = -2i(4j-i a 2(j+i) _ h -a-)\Pa) yields the 
XX super-Hamiltonian 



H 



xx 



(4j a 2i+i + a\ j _ 1 a 2j + 2 - h.a 



(4) 



This Hamiltonian corresponds to two uncoupled chains of 
free fermions. The total number of adjoint fermions, 



E 



2N 



it 



is conserved. Note that the anisotropy A in 



the original XXZ Hamiltonian would introduce recombination 
and pair creation across the chains. Although the above map- 
ping is non local, operators acting only left of a given site j 
will be mapped to fermions that are again only left of this very 
site. So the OSRE of the original XX-model will be the same 
as the corresponding state vector Renyi entropy of two uncou- 
pled chains of free fermions. Thus we have to calculate the en- 
tanglement dynamics of the two uncoupled chains with an ini- 
tial state given by the operator in questions to get the OSRE in 
the XX model. The key point is that local operators are equiv- 
alent to very special, simple initial states in the corresponding 
fermion chains. We here have to distinguish between finite 
index operators (those that involve only a finite number of ad- 
joint fermions after the mapping) and infinite index operators 
(involving a number proportional to the system size L). An 



example of the first kind is |cr| 



Hj-i 



Exam- 



ples of the second kind arise either from local operators like 



W) 



T2(j-1) 



or non local ones like 



nizi^)=i j - l (uli 

We proceed by showing that the bi-partite Renyi entropy 
52 for a system of free fermions in ID is strictly related to 
the number fluctuations in any one of the two partitions as- 
suming a fixed total number. We can assume that the ini- 
tial state of the fermions corresponding to the local opera- 
tors of interest is a Gaussian state. Due to the free evolu- 
tion it remains Gaussian and can be transformed into a prod- 
uct form q = Qj where the Qj correspond to site j and 

have eigenvalues £ 3 > \Vj\ ^ !■ The square of the vari- 
ance of the total particle number in each partition is then 
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AiVf UJ]. On the other hand 



E^(i-^ 2 )/4 

- log 2 Trg 2 = - £\ log 2 (l - (1 - ?/ 2 )/2). Using 

one 



lh 2^ < - l°g 2 (1 - a) < T^. where < x < ± 



obtains 



— AN 2 > 5 2 > — AN 2 
In 2 ~ In 2 



(5) 



For finite index operators we find saturation as can be seen 
in Fig. Q] This reflects the fact, that there is only a finite num- 
ber M of free particles present in both chains together. Thus a 
finite x °f 2 M yields the exact solution |[l6ll for all times |24|. 
For infinite index operators we observe logarithmic growth of 
the OSRE, see Fig. Q] While the infinite number of involved 
adjoint fermions may suggest a linear growth of AN 2 , this is 
not the case as can be understood in the following way: The 
super state corresponding to a infinite index operator like \F) 
(a finite size example of which is shown in Fig. [T]) is filled up 
completely with fermions in the left part of the chains. Inside 
these regions the Pauli principle prevents hopping of fermions 
and thus only particles at the edge where the effective band- 
insulator is connected to the vacuum can move and fill the 
empty parts of the double chain. For the half filled chain An- 
tal et al. have shown that AN 2 = (In t + D)/2tt 2 in the limit 
of large t with a known constant D > 11711 . Other infi- 
nite index operators that result in a initial occupation of the 
two chains different from that of | F) only on a finite num- 
ber of sites show the same logarithmic long time behavior of 
the OSRE, see \a + ) at a single site, also shown in Fig. [T] 
This explains the dynamics of the OSRE in the XX model as 
a generic example of an integrable model that can be mapped 
to free fermions. 

non-integrable models: We now show that there is an- 
other class of systems and operators which may allow an ef- 
ficient simulation of dynamics in the HP. We construct an 
upper bound for the OSRE S a , a > 1, in terms of the 
infinite-temperature auto-correlation function (ITAC). With- 
out loss of generality we assume a normalized operator, i.e. 



-Tr 



= 1, where d is the local dimension of the chain. 
With respect to a splitting of the chain of length L into two 
parts here and below all A act on the sub chain A of length 
La and all B on B of length Lb- Any operator can be rep- 
resented as O (t) = J2 m ,n (*) Am §5 B n with orthonor- 



mal bases 



-Tr 



An A r 



Bl B„ 



A 



= d„ 

is a matrix and its singular values \/A n , (Ai > A2 > ... 
are the eigenvalues of k) are coefficients of a Schmidt de- 
composition 0(t) = J2n=i AiAjf*) ® S n (t), where now 
the Schmidt rank x is at most d 2min( - LA > LB \ This allows to 
express the infinite-temperature auto-correlation function in 
terms of Schmidt coefficients. We find for a > 1 



6\t)d 



= |Tr[At(t)A(0)]| 
x 

< ^ XkXk (°) (6) 



1_A / X 



vfc=l 



(7) 



(8) 



In © we made use of von Neumann's trace inequality (see 
e.g. IU8I0 . Furthermore Jensen's inequality can be used be- 
cause x° is a concave function in x. Finally I© is true by 
the Cauchy-Schwarz inequality. We thus obtain the follow- 
ing estimate for Renyi entropies, assuming an initial product 
operator, Try/ k(0) = 1, for simplicity: 



Sa < 1 log 2 

1 — a 



6\t)d 



T=oo 



fora>l. (9) 



If the ITAC decays with a power law or slower in time, S a will 
grow at most logarithmically for a > 1. The ITAC has been 
studied over decades in condensed matter physics as it is mea- 
sured in nuclear magnetic resonance and neutron scattering 
experiments in magnetic spin chains. While not proofed rigor- 
ously, it is believed that the Blombergen-de Gennes conjecture 
1 191 of spin diffusion holds: If Ylj=i Oj * s a conserved quan- 



tity, then the ITAC of Oj will show diffusive behavior (i.e. 



~ 1/v t in ID). To our knowledge there is no counter exam- 
ple except for integrable models, where this diffusive behavior 



can turn into a ballistic one (i.e. ~ \/t in ID) 11201 12111 . Nev- 
ertheless it always remains slower than exponential. We con- 
clude that in the HP TEBD we can expect 52 to grow at most 
logarithmically in time, even if the model is non-integrable, 
if the initial operator belongs to a conservation law (for inte- 
grable systems there is an infinite number of those, but one 
is sufficient). This in turn indicates that an efficient classical 
simulation should be possible for large times. 

The spin-1 XXZ chain is an example of a non-integrable 
system, although extension to additional hi ghe r-order nonlin- 



fe=i 



ear terms may turn it into an integrable one 122L 12311 . However 
the total z-magnetization Y^j=i &j i s conserved. This con- 
servation law will lead to a logarithmic scaling of S2 for a z . 
Fig. |2]shows numerical indication for this. It should be noted 
that the spin-1 model is computationally much harder than the 
spin-i model since the local Hilbert space dimension is in- 
creased. Although we do observe logarithmic scaling of the 
OSRE corresponding to a z ' , the prefactor is large, such that 
we can not go to far in time. The plot shows data for differ- 
ent matrix dimension x U P to the point where the cutoff error 
becomes substantial. A clear tendency is visible: On the loga- 
rithmic scale S2 approaches a straight line, while in the linear 
plot a sub-linear scaling is evident. This is consistent with 
the expected logarithmic scaling of the OSRE. For <r + Fig. [2] 
shows logarithmic scaling of S<x only for A = 1 because only 
then the total x- and y-magnetization are also conserved. Oth- 
erwise it indicates linear growth of S2 with time. We can un- 
derstand this now as a direct consequence of the Blombergen- 
de Gennes conjecture, which predicts a power law rather than 
an exponentially decaying ITAC in the isotropic case (see in- 
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FIG. 2: (Color online) OSRE dynamics for the 40 site spin-1 XXZ model for a split in the center. Dotted, dashed, and solid lines indicate 
simulations using \ = 300, \ = 500, and \ = 1000 respectively. The left panel features a logarithmic, the right panel a linear time scale. The 
curves show clear indication of the predicted long-time scaling. The insets show two corresponding ITAC curves for a^o ( n °te the logarithmic 
vertical scaling). 



sets of Fig. |2]). We note, that the regular and chaotic Hamil- 
tonians used in the original work by Prosen and Znidaric J3l 
have also been investigated numerically. The OSRE shows 
the expected time dependence, i.e. logarithmic scaling in the 
regular and linear scaling in the chaotic case. 

From the numerical results we can also extract the von Neu- 
mann entropy as a function of time. It scales exactly as Sa in 
the spin-i model for all operators we looked at. The results 
are not conclusive in the spin-1 case however, since the de- 
pendence on the matrix dimension \ use d i n the simulations 
is much stronger. At least they do not contradict the presump- 
tion, that again the scaling is the same as for 5*2. 

In summary we have given a simple explanation of the at 
most logarithmic time dependence of the OSRE 52 for the 
spin- 1/2 XX model as a generic integrable model that can 
be mapped to free fermions. The operator dynamics in that 
model is equivalent to two uncoupled chains of free fermions 
with an initial state corresponding to the operator under con- 
sideration. For local operators these initial states are rather 
simple. E.g. an operator <rj corresponds to a single fermion 
in each chain. We have shown that the bi-partite OSRE 52 
is strictly related to the fluctuations of the fermion number in 
the two partitions, which in turn allowed a simple understand- 
ing of the entropy dynamics. We have shown furthermore that 
for any model, integrable or not, 52 in OS can be bound by 
the infinite-temperature auto-correlation function of the con- 
sidered operator. This in turn means that for systems and ob- 
servables for which the Blombergen-de Gennes conjecture of 
spin diffusion holds, an at most logarithmic growth of the OS 
entanglement is expected. The latter applies e.g. for local 
operators that constitute a global conservation law. 
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cial support, and the Erwin Schrodinger Institute, Vienna, for 
hospitality. 
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